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ABSTRACT 

We investigate the question of whether ambipolar diffusion (ion-neutral drift) 
determines the smallest length and mass scale on which structure forms in a 
turbulent molecular cloud. We simulate magnetized turbulence in a mostly neu- 
tral, uniformly driven, turbulent medium, using a three-dimensional, two-fluid, 
magnetohydrodynamics (MHD) code modified from Zeus-MP. We find that sub- 
stantial structure persists below the ambipolar diffusion scale because of the 
propagation of compressive slow MHD waves at smaller scales. Contrary to sim- 
ple scaling arguments, ambipolar diffusion thus does not suppress structure below 
its characteristic dissipation scale as would be expected for a classical diffusive 
process. We have found this to be true for the magnetic energy, velocity, and 
density. Correspondingly, ambipolar diffusion leaves the clump mass spectrum 
unchanged. Ambipolar diffusion appears unable to set a characteristic scale for 
gravitational collapse and star formation in turbulent molecular clouds. 

Subject headings: ISM: clouds - ISM: kinematics and dynamics - ISM: magnetic 
fields - MHD - stars: formation - turbulence 



1. Introduction 



Molecular clouds are turbulent, with linewidths indicating highly supersonic motions 
(Zuckerman & Palmer 1974), and magnetized, with magnetic energies in or near equipartition 
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with thermal energy (Crutcher 1999). They have low ionization fractions (Elmegreen 1979) 
leading to imperfect coupling of the magnetic field with the gas. Molecular clouds are the sites 
of all known star formation, so characterizing the properties of this non-ideal, magnetized 
turbulence appears central to formulating a theory of star formation. 

The drift of an ionized, magnetized gas through a neutral gas coupled to it by ion-neutral 
collisions is known by astronomers as ambipolar diffusion (AD) and by plasma physicists as 
ion-neutral drift. It was first proposed in an astrophysical context by Mestel & Spitzer (1956) 
as a mechanism for removing magnetic flux and hence magnetic pressure from collapsing 
protostellar cores in the then-novel magnetic field of the Galaxy. However, more recently, as 
turbulence has regained importance in the theory of star formation, AD has been invoked as 
a source of dissipation for magnetic energy in the turbulent magnetohydrodynamic (MHD) 
cascade and thus a characteristic length scale for the star formation process (e.g., Tassis 
& Mouschovias 2004). This is due to its well-known ability to damp certain families of 
linear MHD waves (Kulsrud & Pearce 1969; Ferriere, Zweibel, & Shull 1988; Balsara 1996). 
However, as Balsara (1996) pointed out, AD does allow slow modes to propagate undamped. 

A brief calculation suggests that AD should be the most important dissipation mecha- 
nism in molecular clouds. AD can be expressed as an additional force term in the momentum 
equation for the ions 

F in = PipnlAD(Vn ~ Vi), (1) 

and an equal and opposite force F n i = —Fi n in the neutral momentum equation, where 
and p n are the ion and neutral densities and ^ad — 9.2 x 10 13 cm 3 s _1 g _1 is the collisional 
coupling constant (Draine, Roberge, & Dalgarno 1983; Smith & Mac Low 1997). 

The effect of ion-neutral drift on the magnetic field can be simply expressed in the 
strong coupling approximation (Shu 1983) that neglects the momentum and pressure of the 
ion fluid and equates the collisional drag force on the ions F in with the Lorentz force, 

, (V x B) x B 

-piPnlAD^i - V n ) = . (2) 

47T 

Brandenburg & Zweibel (1994) note that by substituting equation (2) into the induction 
equation for the ions, one arrives at 



d t B = V x 

where 



(v n x B) + (. VxB )- B B - (rj + VAD )V x B 

^PiPnlAD 



(3) 



Vad = -. (4) 

^PiPnlAD 

is the ambipolar diffusivity and r\ is the Ohmic diffusivity. However, dissipation is not the 
only contribution of AD to the induction equation. Given that AD tends to force magnetic 
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fields into force-free states (Brandenburg et al. 1995; Zweibel & Brandenburg 1997) with 
(V x B) x B = 0, it should come as little surprise that the (V x B) • B term must be given 
proper consideration. 

We can approximate the scale ids below which dissipation dominates turbulent structure 
for a given diffusivity r\ in at least two ways. The first is commonly used in the turbulence 
community. It is to equate the driving timescale 

Tdr = L dr /V dr , (5) 

where L dr is the driving wavelength and v dr is the rms velocity at that wavelength, with the 
dissipation timescale T ds = £ ds /f], and solve for £ ds . The second method was suggested by 
Balsara (1996) and Zweibel & Brandenburg (1997) and advocated by Klessen et al. (2000). 
It is to estimate the length scale at which the Reynolds number associated with a given 
dissipation mechanism becomes unity. The Reynolds number for ion-neutral drift can be 
defined as 

Rad = , (6) 

Vad 

where V is a characteristic velocity. This method requires setting Rad to one and solving 
for L = ids to find 

£ad = , ^ „ ■ (7) 

Klessen et al. (2000) show that by adopting values characteristic of dense molecular 
clouds, a magnetic field strength B = 10i?io /iG, ionization fraction x = 10~ 6 X6, neutral 
number density n n = 10 3 n 3 cm -3 , mean mass per particle \i = 2.36m H where m# is the hy- 
drogen mass, such that p n = /in n , and the above value for the ion- neutral coupling constant, 
the length scale at which AD is important is given by 

l AD = (0.04pc) Bl ° (8) 

MaXqU 3 

where Ma = V/va is the Alfven Mach number. By contrast, Ohmic dissipation acts only at 
far smaller scales, £ v ~ 10~ 13 pc (Zweibel & Brandenburg 1997). 

For our purposes, we use the Reynolds number method and choose V = vrms, the RMS 
velocity. Although we use Reynolds numbers, we find that using the timescale method has 
no effect on our results. 

Previous three-dimensional numerical studies of turbulent ion-neutral drift have used 
the strong coupling approximation (Padoan, Zweibel, & Nordlund 2000). This by definition 
renders simulations unable to reach below Rad ~ 1, and thus into the dissipation region. 
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In this paper, we present runs in which we vary the ambipolar diffusion coupling con- 
stant, and thus £ad- We find a surprising lack of dependence of the spectral properties on the 
strength of the ambipolar diffusivity. In particular, no new dissipation range is introduced 
into the density, velocity or magnetic field spectra by ambipolar diffusion, nor is the clump 
mass spectrum materially changed. 

2. Numerical Method 

We solve the two-fluid equations of MHD using the ZEUS-MP code (Norman 2000) 
modified to include a semi-implicit treatment of ion- neutral drift. ZEUS-MP is the domain- 
decomposed, parallel version of the well-known shared memory code ZEUS-3D (Clarke & 
Norman 1994). Both codes follow the algorithms of ZEUS-2D (Stone & Norman 1992a,b), 
including van Leer (1977) advection, and the constrained transport method of characteristics 
(Evans & Hawley 1988; Hawley & Stone 1995) for the magnetic fields. We add an additional 
neutral fluid and collisional coupling terms to both momentum equations. Because ion- 
neutral collisions constitute a stiff term, we evaluate the momentum equations using the 
semi-implicit algorithm of Mac Low & Smith (1997). We also include an explicit treatment 
of Ohmic diffusion by operator splitting the induction equation (Fleming, Stone, & Hawley 
2000). 

We ignore ionization and recombination, assuming that such processes take place on 
timescales much longer than the ones we are concerned with. This means that ions and 
neutrals are separately conserved. Furthermore, we assume that both fluids are isothermal 
and at the same temperature, thus sharing a common sound speed c s . 

2.1. Initial Conditions and Parameters 

All of our runs are on three-dimensional Cartesian grids with periodic boundary condi- 
tions in all directions. 

The turbulence is driven by the method detailed in Mac Low (1999). Briefly, we generate 
a top hat function in Fourier space between 1 < \k\ < 2. The amplitudes and phases of each 
mode are chosen at random, and once returned to physical space, the resulting velocities 
are normalized to produce the desired RMS velocity, unity in our case. At each timestep, 
the same pattern of velocity perturbations is renormalized to drive the box with a constant 
energy input (E — 1.0 for all simulations) and applied to the neutral gas. 

Our isothermal sound speed is c s = 0.1, corresponding to an initial RMS Mach number 
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M — 10. The initial neutral density p n is everywhere constant and set to unity. The magnetic 
field strength is set by requiring that the initial ratio of gas pressure to magnetic pressure 
be everywhere (3 = 8irc 2 s p/B 2 = 0.1; its direction lies along the z-axis. 

Although our semi-implicit method means that the timestep is not restricted by the 
standard Courant condition for diffusive processes (that is, oc [Aa;] 2 ), the two-fluid model is 
limited by the Alfven timestep for the ions. This places strong constraints on the ionization 
fraction (x = rii/n n ) we can reasonably compute. We therefore adapt a fixed fraction of 
x — 0.1 for our simulations. While this fraction is certainly considerably higher than the 
10~ 4 -10~ 9 typical of molecular clouds, the ionization fraction only enters the calculation in 
concert with the collisional coupling constant 7 ad- Thus, we are able to compensate for the 
unrealistically high ionization fraction by adjusting ^ad accordingly. 

We present four runs, two with AD, one with Ohmic diffusion, and one ideal MHD run 
(see Table 1). For the AD runs, we vary the collisional coupling constant in order to change 
the diffusivity. 

Our results are reported for a resolution of 256 3 at time t = 0.125t s = 2.5 where 
t s = 20 is the sound crossing time for the box. This exceeds by at least 30% the turbulent 
crossing time over the driving scale computed from equation (5), and tabulated in Table 1. 
Our computation of r^r is done for = 1, the maximum driving wavelength. Vestuto, 
Ostriker, & Stone (2003) note that r dr is the relevant timescale for the formation of nonlinear 
structures. Furthermore, we find from studies performed at 128 3 out to t — 0.3t s that 0.125t s 
is enough time to reach a steady state in energy. 

3. Results 

Figure 1 shows cuts of density perpendicular and parallel to the mean magnetic field. 
For the ambipolar diffusion runs, we show the total density p = Pi + p n - The morphology of 
density enhancements in the different runs appears similar, giving a qualitative suggestion 
of the quantitative results on clump mass spectra discussed next. 

3.1. Clump mass spectrum 

We wish to understand whether AD determines the smallest scale at which clumps 
can form in turbulent molecular clouds. Determining structure within molecular clouds has 
proved difficult in both theory and observation. Molecular line maps (eg, Falgarone et al 
1992) show that for all resolvable scales, the density fields of clouds is made up of a hierarchy 
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of clumps. Furthermore, the identification of clumps projected on the sky with physical 
volumetric objects is questionable (Ostriker, Stone, & Gammie 2001; Ballesteros-Paredes & 
Mac Low 2002). 

Nonetheless, density enhancements in a turbulent flow likely provide the initial condi- 
tions for star formation. To clarify the effects of different turbulent dissipation mechanisms 
on the clump mass spectrum, we study our three dimensional simulations of turbulence with- 
out gravity. By using the clumpfind algorithm (Williams, de Geus, & Blitz 1994) on the 
density field to identify contiguous regions of enhanced density, we can construct a clump 
mass spectrum (Fig. 2). Although such methods are parameter-sensitive when attempting to 
draw comparisons to observed estimates for the clump-mass spectrum (Ballesteros-Paredes 
& Mac Low 2002), we are only interested in using the mass spectrum as a point of comparison 
between runs with different dissipative properties. 

For this section, we dimensionalize our density field following Mac Low (1999), with a 
length scale L' = 0.5 pc, and mean density scale p' = 10 4 (2m#) g cm -3 in order to present 
results in physical units relevant to star formation. 

We search for clumps above a density threshold set at 5(p) (where in the AD cases 
p = Pi + p n ) and bin the results by mass to produce a clump-mass spectrum. Figure 2 shows 
that while Ohmic diffusion has a dramatic effect on the number of low-mass clumps, AD 
has nearly none. Although there are small fluctuations around the hydrodynamic spectrum, 
there is no systematic trend with increasing strength of AD. This result suggests that AD 
does not control the minimum mass of clumps formed in turbulent molecular clouds. 

3.2. Magnetic Energy and Density Spectra 

The lack of an effect on the clump mass spectrum can be better understood by examining 
the distribution of magnetic field and density. 

AD produces no evident dissipation range in the magnetic energy spectrum. As seen 
in Figure 3, for two different values of ambipolar diffusivity tjad, the power spectrum of 
magnetic field retains the shape of the ideal run. For comparison, we have also plotted the 
run with Ohmic diffusion. While the expected dissipation wavenumbers (determined in both 
cases by the Reynolds number method mentioned above) of the t]ad = 0.275 and rj = 0.250 
runs are very similar, the effect of Ohmic diffusion is quite apparent in the declining slope 
of the magnetic energy spectrum, in contrast to AD. 

The total power does decrease as the ambipolar diffusivity tjad increases. Because we 
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drive only the neutrals, this could be interpreted as magnetic energy being lost during the 
transfer of driving energy from the ions to the neutrals via the coupling. However, we 
performed a simulation in which both ions and neutrals were driven with the same driving 
pattern and found almost no difference in the power spectra from our standard (neutral 
driving only) case. 

We instead suspect that the decline in total magnetic energy occurs because AD does 
damp some families of MHD waves, notably Alfven waves (Kulsrud & Pearce 1969), even 
though it does not introduce a characteristic damping scale. 

In order to demonstrate this, the flow will need to be decomposed into its constituent 
MHD wave motions at each point in space. Such a technique has been used before by Maron 
& Goldreich (2001) for incompressible MHD turbulence and by Cho & Lazarian (2002) for 
compressible MHD turbulence. The technique used by Cho & Lazarian (2002) decomposes 
wave motions along a mean field assumed to be present. However, because the local field 
is distorted by the turbulence and thus not necessarily parallel to the mean, a mean-field 
decomposition tends to spuriously mix Alfven and slow modes (Maron & Goldreich 2001). If 
the local field line distortion is great enough, the decomposition must be made with respect to 
the local field, a much more demanding proceedure. Although wave decomposition analysis 
is outside the scope of this paper, it remains a fruitful avenue for future research. 

In order to ensure that the lack of spectral features seen in the magnetic spectrum 
(and similarly in the density spectrum) is not an artifact of the limited inertial range in 
our simulations, we ran our tjad = 0.275 (medium collision strength) case at resolutions 
of 64 3 , 128 3 , and 256 3 . Figure 4 demonstrates that increasing the resolution increases the 
inertial range, but does not resolve any noticeable transition to dissipation at the AD length, 
suggesting that our results are not sensitive to the resolution. 

Figure 5 shows the spectrum of the density for all runs. In the case of the AD runs, we 
use the sum of the neutral and ion density. 

The density spectrum peaks at small scale in compressible turbulence (Joung & Mac 
Low 2005). Varying the ambipolar diffusivity by a factor of two makes little systematic 
difference to the shape of the density spectrum. It seems clear that although there are only 
slight differences in the density spectrum due to varying magnetic diffusivities, the density 
spectrum is not a particularly good indicator of underlying clump masses. 

Note that we use for the density spectrum the Fourier transform of the density field 
rather than its square, which in the case of the magnetic field yields the one-point correlation 
function (or power spectrum) of the magnetic energy. 
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4. Discussion 

Supersonic turbulence performs a dual role in its simultaneous ability to globally sup- 
port a molecular cloud against gravity while at the same time producing smaller density 
enhancements that can sometimes gravitationally collapse (Klessen et al. 2000). While our 
simulations do not include gravity, it is clear that AD does not set a characteristic scale 
to the density field below which MHD turbulence is unable to further influence structure 
formation. 

One of the main motivations of this study was to verify the claim made by, for exam- 
ple, Klessen et al. (2000) that AD sets the minimum mass for clumps in molecular cloud 
turbulence. However, it appears that AD is unable to set this scale, because of its selective 
action on different MHD waves. We do note that AD can occasionally help form magneto- 
hydrostatic objects in MHD turbulence, but this is not a dominant pathway, as shown by 
Vazquez-Semadeni et al. (2005). Although Ohmic diffusion has little trouble inhibiting low 
mass clump formation, it never reaches significant values at the densities where molecular 
clumps form. 

This opens up other possibilities for the physical mechanisms determining the smallest 
scale fluctuations occurring in molecular clouds. An attractive option is the sonic-scale 
argument of Vazquez-Semadeni, Ballesteros-Paredes, & Klessen (2003), in which the length 
scale at which turbulent motions become incompressible, with Mach numbers dropping well 
below unity, determines where turbulence ceases to have an effect on the pre-stellar core 
distribution, and thus determines the minimum mass scale. 
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Fig. 1. — Random cuts of density p parallel and perpendicular to the magnetic field for each 
of three runs of varying ambipolar diffusivity tjad ■ Each image is scaled to its own minimum 
and maximum, enhancing structural features. For AD runs, p = Pi + p n 
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Fig. 2. — Clump mass spectrum measuring the number of clumps of a given mass for one ideal 
MHD run (labeled tjad = 0), two AD runs (tjad — 0.275,0.525) and one Ohmic dissipation 
run (rj = 0.250). Compare the lack of effect of AD to the significant decrease in the number 
of low mass clumps for the Ohmic diffusion case. 
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Fig. 3. — Magnetic energy spectra for the same runs as Figure 2. The vertical lines represent 
the wavenumber at which the AD or Ohmic Reynolds number crosses unity. 
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Fig. 4. — Magnetic energy spectra for three runs of varying resolution from 64 3 to 256 3 . 
Increased resolution shows no effects at the AD wavenumber given by the vertical line. 
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Fig. 5. — Fourier transform of density field for the same runs as Figure 2. Note that for the 
AD runs, this plot shows the total density, p = Pi + p n - 
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Table 1. Models 



Run 


diffusivity 


IAD 


Vad 


V 






Tdr 


Al 


AD 


8 


0.275 





0.603 


0.526 


1.66 


A2 


AD 


4 


0.575 





0.615 


0.501 


1.63 





Ohmic 








0.250 




0.577 


1.73 


I 















0.630 


1.59 



